{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7a1cfa3d-f83d-4c10-a7ab-fb986c297dfb",
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import statsmodels.formula.api as smf\n",
    "import matplotlib.pyplot as plt\n",
    "from matplotlib.lines import Line2D\n",
    "from IPython.display import display, Markdown\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "92a9db3b-af88-4a37-adcf-6bdcc59c46f7",
   "metadata": {},
   "outputs": [],
   "source": [
    "# =========================\n",
    "# 1. PATHS\n",
    "# =========================\n",
    "DATA_DIR = r\"E:\\\"\n",
    "DATA_FILE = os.path.join(DATA_DIR, \"df_1.xlsx\")\n",
    "SHEET_NAME = \"Data\"\n",
    "\n",
    "OUT_DIR = os.path.join(DATA_DIR, \"outputs_publication\")\n",
    "FIG_DIR = os.path.join(OUT_DIR, \"figures\")\n",
    "TAB_DIR = os.path.join(OUT_DIR, \"tables\")\n",
    "RES_DIR = os.path.join(OUT_DIR, \"full_results\")\n",
    "\n",
    "for p in [OUT_DIR, FIG_DIR, TAB_DIR, RES_DIR]:\n",
    "    os.makedirs(p, exist_ok=True)\n",
    "\n",
    "# =========================\n",
    "# 2. LOAD DATA\n",
    "# =========================\n",
    "df = pd.read_excel(DATA_FILE, sheet_name=SHEET_NAME).copy()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "bca2d2c9-c049-42de-b1ec-b676aa42d653",
   "metadata": {},
   "outputs": [],
   "source": [
    "# =========================\n",
    "# 3. VARIABLE SETUP\n",
    "# =========================\n",
    "controls = [\n",
    "    \"Q26_age\",\n",
    "    \"Q27_gender\",\n",
    "    \"Q28_region\",\n",
    "    \"Q29_education\",\n",
    "    \"Q30_income\",\n",
    "    \"Q31_marital\",\n",
    "    \"Q32_union\",\n",
    "    \"Q33_employment\",\n",
    "    \"Q34_homeownership\",\n",
    "    \"Q38_subjective_class\",\n",
    "    \"Q39_ideology\",\n",
    "]\n",
    "\n",
    "emotion_items = [\"Q6\", \"Q7\", \"Q8\", \"Q9\", \"Q10\"]\n",
    "attribution_items = [\"Q11\", \"Q12\", \"Q13\", \"Q14\"]\n",
    "policy_items = [\"Q15\", \"Q16\", \"Q17\", \"Q18\", \"Q19\"]\n",
    "institution_items = [\"Q20\", \"Q21\", \"Q22\", \"Q23\"]\n",
    "\n",
    "emotion_index = \"emotional_burden_index\"\n",
    "attribution_index = \"attribution_index\"\n",
    "policy_index = \"policy_index\"\n",
    "institution_index = \"institutional_index\"\n",
    "asset_index = \"asset_index\"\n",
    "\n",
    "categorical_controls = [\n",
    "    \"Q27_gender\",\n",
    "    \"Q28_region\",\n",
    "    \"Q29_education\",\n",
    "    \"Q30_income\",\n",
    "    \"Q31_marital\",\n",
    "    \"Q32_union\",\n",
    "    \"Q33_employment\",\n",
    "    \"Q34_homeownership\",\n",
    "]\n",
    "\n",
    "continuous_controls = [\n",
    "    \"Q26_age\",\n",
    "    \"Q38_subjective_class\",\n",
    "    \"Q39_ideology\",\n",
    "]\n",
    "\n",
    "item_labels = {\n",
    "    \"Q6\":  \"DV: Stress\",\n",
    "    \"Q7\":  \"DV: Anxiety\",\n",
    "    \"Q8\":  \"DV: Exhaustion\",\n",
    "    \"Q9\":  \"DV: Frustration\",\n",
    "    \"Q10\": \"DV: Pressure\",\n",
    "    \"Q11\": \"DV: Social conditions\",\n",
    "    \"Q12\": \"DV: Difficult individually\",\n",
    "    \"Q13\": \"DV: Structural roots\",\n",
    "    \"Q14\": \"DV: Not purely personal\",\n",
    "    \"Q15\": \"DV: Gov. responsibility\",\n",
    "    \"Q16\": \"DV: Support expansion\",\n",
    "    \"Q17\": \"DV: Protection\",\n",
    "    \"Q18\": \"DV: Spending support\",\n",
    "    \"Q19\": \"DV: Public matter\",\n",
    "    \"Q20\": \"DV: Trust decline\",\n",
    "    \"Q21\": \"DV: Institutional inadequacy\",\n",
    "    \"Q22\": \"DV: Perceived unfairness\",\n",
    "    \"Q23\": \"DV: Political resolution\",\n",
    "    emotion_index: \"DV: Emotional burden index\",\n",
    "    attribution_index: \"DV: Attribution index\",\n",
    "    policy_index: \"DV: Policy index\",\n",
    "    institution_index: \"DV: Institutional index\",\n",
    "}\n",
    "\n",
    "# =========================\n",
    "# 4. STYLE\n",
    "# =========================\n",
    "plt.rcParams.update({\n",
    "    \"font.family\": \"DejaVu Sans\",\n",
    "    \"font.size\": 12,\n",
    "    \"axes.labelsize\": 12,\n",
    "    \"axes.titlesize\": 16,\n",
    "    \"xtick.labelsize\": 12.5,\n",
    "    \"ytick.labelsize\": 12.5,\n",
    "    \"legend.fontsize\": 9,   # reduced\n",
    "    \"figure.dpi\": 150,\n",
    "    \"savefig.dpi\": 600,\n",
    "    \"axes.spines.top\": False,\n",
    "    \"axes.spines.right\": False,\n",
    "})\n",
    "\n",
    "MARKER_SIZE = 7.5\n",
    "LINE_WIDTH = 1.7\n",
    "TITLE_PAD = 14\n",
    "ANNOT_SIZE = 9.2\n",
    "LEGEND_HANDLE_LEN = 1.7\n",
    "\n",
    "legend_handles = [\n",
    "    Line2D([0], [0], marker='o', linestyle='None', markersize=5.5,\n",
    "           label='Estimate', color='black'),\n",
    "    Line2D([0], [0], linestyle='-', linewidth=1.3, color='black',\n",
    "           label='95% CI')\n",
    "]\n",
    "\n",
    "# =========================\n",
    "# 5. HELPERS\n",
    "# =========================\n",
    "def build_formula(depvar, include_emotion=False, include_asset_interaction=False):\n",
    "    rhs = [\"other\"]\n",
    "\n",
    "    if include_emotion:\n",
    "        rhs.append(emotion_index)\n",
    "\n",
    "    if include_asset_interaction:\n",
    "        rhs += [asset_index, f\"other:{asset_index}\"]\n",
    "\n",
    "    for var in continuous_controls:\n",
    "        rhs.append(var)\n",
    "\n",
    "    for var in categorical_controls:\n",
    "        rhs.append(f\"C({var})\")\n",
    "\n",
    "    return f\"{depvar} ~ \" + \" + \".join(rhs)\n",
    "\n",
    "\n",
    "def run_model(depvar, include_emotion=False, include_asset_interaction=False):\n",
    "    formula = build_formula(\n",
    "        depvar=depvar,\n",
    "        include_emotion=include_emotion,\n",
    "        include_asset_interaction=include_asset_interaction\n",
    "    )\n",
    "    model = smf.ols(formula, data=df).fit(cov_type=\"HC3\")\n",
    "    return model\n",
    "\n",
    "\n",
    "def extract_term(model, term, depvar, dep_label):\n",
    "    coef = model.params.get(term, np.nan)\n",
    "    se = model.bse.get(term, np.nan)\n",
    "    ci_low, ci_high = model.conf_int().loc[term] if term in model.params.index else (np.nan, np.nan)\n",
    "    pval = model.pvalues.get(term, np.nan)\n",
    "    return {\n",
    "        \"depvar\": depvar,\n",
    "        \"dep_label\": dep_label,\n",
    "        \"term\": term,\n",
    "        \"coef\": coef,\n",
    "        \"se\": se,\n",
    "        \"ci_low\": ci_low,\n",
    "        \"ci_high\": ci_high,\n",
    "        \"p_value\": pval,\n",
    "        \"n\": int(model.nobs),\n",
    "        \"r2\": model.rsquared,\n",
    "        \"adj_r2\": model.rsquared_adj,\n",
    "    }\n",
    "\n",
    "\n",
    "def style_axis(ax, xlab=None, ylab=None):\n",
    "    ax.axvline(0, color=\"grey\", linestyle=\"--\", linewidth=1.0, alpha=0.8, zorder=0)\n",
    "    ax.tick_params(axis='both', which='major', length=5, width=1.0, direction='out')\n",
    "    ax.tick_params(axis='both', which='minor', length=3, width=0.8, direction='out')\n",
    "    ax.spines[\"left\"].set_linewidth(1.0)\n",
    "    ax.spines[\"bottom\"].set_linewidth(1.0)\n",
    "    if xlab is not None:\n",
    "        ax.set_xlabel(xlab)\n",
    "    if ylab is not None:\n",
    "        ax.set_ylabel(ylab)\n",
    "\n",
    "\n",
    "def coef_plot(ax, results_df, subtitle, xlim=None, annotate_coef=True):\n",
    "    plot_df = results_df.copy().reset_index(drop=True)\n",
    "    plot_df = plot_df.iloc[::-1].reset_index(drop=True)\n",
    "    y = np.arange(len(plot_df))\n",
    "\n",
    "    ax.hlines(\n",
    "        y=y,\n",
    "        xmin=plot_df[\"ci_low\"],\n",
    "        xmax=plot_df[\"ci_high\"],\n",
    "        color=\"black\",\n",
    "        linewidth=LINE_WIDTH,\n",
    "        zorder=2\n",
    "    )\n",
    "    ax.plot(\n",
    "        plot_df[\"coef\"], y,\n",
    "        \"o\",\n",
    "        color=\"black\",\n",
    "        markersize=MARKER_SIZE,\n",
    "        zorder=3\n",
    "    )\n",
    "\n",
    "    if annotate_coef:\n",
    "        xspan = (xlim[1] - xlim[0]) if xlim is not None else (plot_df[\"ci_high\"].max() - plot_df[\"ci_low\"].min())\n",
    "        offset = 0.08 * xspan\n",
    "        for _, row in plot_df.iterrows():\n",
    "            x = row[\"coef\"]\n",
    "            yy = y[plot_df.index[plot_df[\"dep_label\"] == row[\"dep_label\"]][0]]\n",
    "            label = f\"{row['coef']:.3f}\"\n",
    "            if np.isfinite(x):\n",
    "                if x < 0:\n",
    "                    ax.text(x - offset, yy, label, ha=\"right\", va=\"center\", fontsize=ANNOT_SIZE)\n",
    "                else:\n",
    "                    ax.text(x + offset, yy, label, ha=\"left\", va=\"center\", fontsize=ANNOT_SIZE)\n",
    "\n",
    "    ax.set_yticks(y)\n",
    "    ax.set_yticklabels(plot_df[\"dep_label\"])\n",
    "    ax.set_title(subtitle, loc=\"center\", pad=TITLE_PAD, fontsize=16)\n",
    "    if xlim is not None:\n",
    "        ax.set_xlim(xlim)\n",
    "    style_axis(ax, xlab=\"Coefficient on Other (Self=0, Other=1)\")\n",
    "    ax.legend(\n",
    "        handles=legend_handles,\n",
    "        loc=\"upper left\",\n",
    "        frameon=False,\n",
    "        handlelength=LEGEND_HANDLE_LEN,\n",
    "        borderpad=0.15,\n",
    "        labelspacing=0.3,\n",
    "        handletextpad=0.45\n",
    "    )\n",
    "\n",
    "\n",
    "def save_full_results(models_dict, outfile_excel, outfile_csv_prefix):\n",
    "    rows = []\n",
    "    for name, model in models_dict.items():\n",
    "        conf = model.conf_int()\n",
    "        for term in model.params.index:\n",
    "            rows.append({\n",
    "                \"model\": name,\n",
    "                \"term\": term,\n",
    "                \"coef\": model.params[term],\n",
    "                \"se\": model.bse[term],\n",
    "                \"ci_low\": conf.loc[term, 0],\n",
    "                \"ci_high\": conf.loc[term, 1],\n",
    "                \"p_value\": model.pvalues[term],\n",
    "                \"n\": int(model.nobs),\n",
    "                \"r2\": model.rsquared,\n",
    "                \"adj_r2\": model.rsquared_adj,\n",
    "            })\n",
    "\n",
    "    res_df = pd.DataFrame(rows)\n",
    "    res_df.to_excel(outfile_excel, index=False)\n",
    "    res_df.to_csv(f\"{outfile_csv_prefix}.csv\", index=False, encoding=\"utf-8-sig\")\n",
    "    return res_df\n",
    "\n",
    "\n",
    "def annotate_endpoints(ax, x, y, condition=\"Self\", label_prefix=\"\", fontsize=9.2, y_offset=0.06):\n",
    "    x = np.asarray(x)\n",
    "    y = np.asarray(y)\n",
    "\n",
    "    # Self는 숫자를 점 위에, Other는 점 아래에\n",
    "    if condition == \"Self\":\n",
    "        va_start = \"bottom\"\n",
    "        va_end = \"bottom\"\n",
    "        y_start = y[0] + y_offset\n",
    "        y_end = y[-1] + y_offset\n",
    "    else:\n",
    "        va_start = \"top\"\n",
    "        va_end = \"top\"\n",
    "        y_start = y[0] - y_offset\n",
    "        y_end = y[-1] - y_offset\n",
    "\n",
    "    # first point\n",
    "    ax.text(\n",
    "        x[0], y_start,\n",
    "        f\"{label_prefix}{y[0]:.2f}\",\n",
    "        ha=\"center\", va=va_start,\n",
    "        fontsize=fontsize\n",
    "    )\n",
    "\n",
    "    # last point\n",
    "    ax.text(\n",
    "        x[-1], y_end,\n",
    "        f\"{label_prefix}{y[-1]:.2f}\",\n",
    "        ha=\"center\", va=va_end,\n",
    "        fontsize=fontsize\n",
    "    )\n"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "57c893aa-2c93-4e8a-bc91-bc769c7eb54e",
   "metadata": {},
   "source": [
    "## Figure 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3efce65a-a94b-4fc4-9934-83dd6434ea9f",
   "metadata": {},
   "outputs": [],
   "source": [
    "# =========================\n",
    "# 6. FIGURE 1\n",
    "# =========================\n",
    "fig1_models = {}\n",
    "fig1_rows = []\n",
    "\n",
    "for q in emotion_items:\n",
    "    m = run_model(q)\n",
    "    fig1_models[f\"Figure1_{q}\"] = m\n",
    "    fig1_rows.append(extract_term(m, \"other\", q, item_labels[q]))\n",
    "\n",
    "fig1_panelA = pd.DataFrame(fig1_rows)\n",
    "\n",
    "m_em_idx = run_model(emotion_index)\n",
    "fig1_models[\"Figure1_emotional_index\"] = m_em_idx\n",
    "fig1_panelB = pd.DataFrame([\n",
    "    extract_term(m_em_idx, \"other\", emotion_index, item_labels[emotion_index])\n",
    "])\n",
    "\n",
    "save_full_results(\n",
    "    fig1_models,\n",
    "    outfile_excel=os.path.join(RES_DIR, \"figure1_full_results.xlsx\"),\n",
    "    outfile_csv_prefix=os.path.join(RES_DIR, \"figure1_full_results\")\n",
    ")\n",
    "\n",
    "fig, axes = plt.subplots(\n",
    "    nrows=1, ncols=2,\n",
    "    figsize=(13.0, 6.8),\n",
    "    gridspec_kw={\"width_ratios\": [1.5, 0.95]}\n",
    ")\n",
    "\n",
    "coef_plot(\n",
    "    axes[0],\n",
    "    fig1_panelA,\n",
    "    subtitle=\"(a) Emotional burden items\",\n",
    "    xlim=(-1.1, 0.25),\n",
    "    annotate_coef=True\n",
    ")\n",
    "\n",
    "coef_plot(\n",
    "    axes[1],\n",
    "    fig1_panelB,\n",
    "    subtitle=\"(b) Emotional burden index\",\n",
    "    xlim=(-1.1, 0.25),\n",
    "    annotate_coef=True\n",
    ")\n",
    "\n",
    "plt.tight_layout()\n",
    "fig.savefig(os.path.join(FIG_DIR, \"Figure1_emotional_burden.png\"), bbox_inches=\"tight\")\n",
    "fig.savefig(os.path.join(FIG_DIR, \"Figure1_emotional_burden.pdf\"), bbox_inches=\"tight\")\n",
    "plt.close(fig)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c4b88682-f0f7-407a-aa27-b80991786acc",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "8c768298-f8c1-4abc-977f-3f622bb34402",
   "metadata": {},
   "source": [
    "## Figure 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7b3ba1f1-3e29-4c2b-91ee-ad77468f3b91",
   "metadata": {},
   "outputs": [],
   "source": [
    "def coef_plot(ax, results_df, subtitle, xlim=None, annotate_coef=True):\n",
    "    plot_df = results_df.copy().reset_index(drop=True)\n",
    "    plot_df = plot_df.iloc[::-1].reset_index(drop=True)\n",
    "    y = np.arange(len(plot_df))\n",
    "\n",
    "    ax.hlines(\n",
    "        y=y,\n",
    "        xmin=plot_df[\"ci_low\"],\n",
    "        xmax=plot_df[\"ci_high\"],\n",
    "        color=\"black\",\n",
    "        linewidth=LINE_WIDTH,\n",
    "        zorder=2\n",
    "    )\n",
    "    ax.plot(\n",
    "        plot_df[\"coef\"], y,\n",
    "        \"o\",\n",
    "        color=\"black\",\n",
    "        markersize=MARKER_SIZE,\n",
    "        zorder=3\n",
    "    )\n",
    "\n",
    "    if annotate_coef:\n",
    "        xspan = (xlim[1] - xlim[0]) if xlim is not None else (plot_df[\"ci_high\"].max() - plot_df[\"ci_low\"].min())\n",
    "        offset = 0.08 * xspan\n",
    "        for _, row in plot_df.iterrows():\n",
    "            x = row[\"coef\"]\n",
    "            yy = y[plot_df.index[plot_df[\"dep_label\"] == row[\"dep_label\"]][0]]\n",
    "            label = f\"{row['coef']:.3f}\"\n",
    "            if np.isfinite(x):\n",
    "                if x < 0:\n",
    "                    ax.text(x - offset, yy, label, ha=\"right\", va=\"center\", fontsize=ANNOT_SIZE)\n",
    "                else:\n",
    "                    ax.text(x + offset, yy, label, ha=\"left\", va=\"center\", fontsize=ANNOT_SIZE)\n",
    "\n",
    "    ax.set_yticks(y)\n",
    "    ax.set_yticklabels(plot_df[\"dep_label\"])\n",
    "    ax.set_title(subtitle, loc=\"center\", pad=TITLE_PAD, fontsize=18)\n",
    "    if xlim is not None:\n",
    "        ax.set_xlim(xlim)\n",
    "    style_axis(ax, xlab=\"Coefficient on Other (Self=0, Other=1)\")\n",
    "    ax.legend(\n",
    "        handles=legend_handles,\n",
    "        loc=\"upper left\",\n",
    "        frameon=False,\n",
    "        handlelength=LEGEND_HANDLE_LEN,\n",
    "        borderpad=0.15,\n",
    "        labelspacing=0.01,\n",
    "        handletextpad=0.45\n",
    "    )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "24f910cf-4c06-4a02-9ae0-2b82ce1efdfd",
   "metadata": {},
   "outputs": [],
   "source": [
    "# =========================\n",
    "# 7. FIGURE 2\n",
    "# =========================\n",
    "fig2_models = {}\n",
    "fig2_rows_A = []\n",
    "fig2_rows_B = []\n",
    "fig2_rows_C = []\n",
    "\n",
    "for q in attribution_items:\n",
    "    m = run_model(q)\n",
    "    fig2_models[f\"Figure2_{q}\"] = m\n",
    "    fig2_rows_A.append(extract_term(m, \"other\", q, item_labels[q]))\n",
    "\n",
    "for q in policy_items:\n",
    "    m = run_model(q)\n",
    "    fig2_models[f\"Figure2_{q}\"] = m\n",
    "    fig2_rows_B.append(extract_term(m, \"other\", q, item_labels[q]))\n",
    "\n",
    "for q in institution_items:\n",
    "    m = run_model(q)\n",
    "    fig2_models[f\"Figure2_{q}\"] = m\n",
    "    fig2_rows_C.append(extract_term(m, \"other\", q, item_labels[q]))\n",
    "\n",
    "fig2_A = pd.DataFrame(fig2_rows_A)\n",
    "fig2_B = pd.DataFrame(fig2_rows_B)\n",
    "fig2_C = pd.DataFrame(fig2_rows_C)\n",
    "\n",
    "save_full_results(\n",
    "    fig2_models,\n",
    "    outfile_excel=os.path.join(RES_DIR, \"figure2_full_results.xlsx\"),\n",
    "    outfile_csv_prefix=os.path.join(RES_DIR, \"figure2_full_results\")\n",
    ")\n",
    "\n",
    "all_ci = pd.concat([fig2_A, fig2_B, fig2_C], axis=0)\n",
    "xmin = min(-1.2, np.floor(all_ci[\"ci_low\"].min() * 10) / 10)\n",
    "xmax = max(0.25, np.ceil(all_ci[\"ci_high\"].max() * 10) / 10)\n",
    "\n",
    "fig, axes = plt.subplots(\n",
    "    nrows=1, ncols=3,\n",
    "    figsize=(16.4, 7.4),\n",
    "    gridspec_kw={\"width_ratios\": [1.0, 1.15, 1.0]}\n",
    ")\n",
    "\n",
    "coef_plot(\n",
    "    axes[0],\n",
    "    fig2_A,\n",
    "    subtitle=\"(a) Attribution\",\n",
    "    xlim=(xmin, xmax),\n",
    "    annotate_coef=True\n",
    ")\n",
    "coef_plot(\n",
    "    axes[1],\n",
    "    fig2_B,\n",
    "    subtitle=\"(b) Policy responsibility\",\n",
    "    xlim=(xmin, xmax),\n",
    "    annotate_coef=True\n",
    ")\n",
    "coef_plot(\n",
    "    axes[2],\n",
    "    fig2_C,\n",
    "    subtitle=\"(c) Institutional evaluation\",\n",
    "    xlim=(xmin, xmax),\n",
    "    annotate_coef=True\n",
    ")\n",
    "\n",
    "plt.tight_layout()\n",
    "fig.savefig(os.path.join(FIG_DIR, \"Figure2_political_translation.png\"), bbox_inches=\"tight\")\n",
    "fig.savefig(os.path.join(FIG_DIR, \"Figure2_political_translation.pdf\"), bbox_inches=\"tight\")\n",
    "plt.close(fig)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1502ad28-103b-4a64-8a2c-2f17148bb85b",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "4071ef79-885c-4a3b-a352-ef86f2f46d2f",
   "metadata": {},
   "source": [
    "## Table 1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a5585550-41a1-40a1-ba40-3845f4785688",
   "metadata": {},
   "outputs": [],
   "source": [
    "# =========================\n",
    "# 8. TABLE 1\n",
    "# =========================\n",
    "m1 = run_model(policy_index, include_emotion=False)\n",
    "m2 = run_model(policy_index, include_emotion=True)\n",
    "\n",
    "table1_rows = []\n",
    "for model_name, model in [(\"Model 1\", m1), (\"Model 2\", m2)]:\n",
    "    conf = model.conf_int()\n",
    "    for term in [\"other\", emotion_index]:\n",
    "        if term in model.params.index:\n",
    "            table1_rows.append({\n",
    "                \"model\": model_name,\n",
    "                \"term\": term,\n",
    "                \"coef\": model.params[term],\n",
    "                \"se_hc3\": model.bse[term],\n",
    "                \"ci_low\": conf.loc[term, 0],\n",
    "                \"ci_high\": conf.loc[term, 1],\n",
    "                \"p_value\": model.pvalues[term],\n",
    "                \"n\": int(model.nobs),\n",
    "                \"r2\": model.rsquared,\n",
    "                \"adj_r2\": model.rsquared_adj,\n",
    "            })\n",
    "\n",
    "table1_df = pd.DataFrame(table1_rows)\n",
    "table1_df.to_excel(os.path.join(TAB_DIR, \"Table1_policy_models.xlsx\"), index=False)\n",
    "table1_df.to_csv(os.path.join(TAB_DIR, \"Table1_policy_models.csv\"), index=False, encoding=\"utf-8-sig\")\n",
    "\n",
    "save_full_results(\n",
    "    {\"Table1_Model1\": m1, \"Table1_Model2\": m2},\n",
    "    outfile_excel=os.path.join(RES_DIR, \"table1_full_results.xlsx\"),\n",
    "    outfile_csv_prefix=os.path.join(RES_DIR, \"table1_full_results\")\n",
    ")\n",
    "\n",
    "with open(os.path.join(TAB_DIR, \"Table1_policy_models.txt\"), \"w\", encoding=\"utf-8\") as f:\n",
    "    f.write(\"TABLE 1. Policy regression models\\n\\n\")\n",
    "    f.write(m1.summary().as_text())\n",
    "    f.write(\"\\n\\n\" + \"=\"*90 + \"\\n\\n\")\n",
    "    f.write(m2.summary().as_text())\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "167c89d3-ce83-4354-b3c8-bb9eaa50f4fa",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "690b8c72-5a9d-4a7e-bd27-280985eab839",
   "metadata": {},
   "outputs": [],
   "source": [
    "import os\n",
    "import pandas as pd\n",
    "from IPython.display import display, Markdown\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 1. OUTPUT PATHS\n",
    "# -------------------------------------------------\n",
    "TAB_DIR = os.path.join(OUT_DIR, \"tables\")\n",
    "RES_DIR = os.path.join(OUT_DIR, \"full_results\")\n",
    "os.makedirs(TAB_DIR, exist_ok=True)\n",
    "os.makedirs(RES_DIR, exist_ok=True)\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 2. DV LIST\n",
    "# -------------------------------------------------\n",
    "table1_dvs = [\n",
    "    attribution_index,\n",
    "    policy_index,\n",
    "    institution_index\n",
    "]\n",
    "\n",
    "dv_pretty = {\n",
    "    attribution_index: \"Attribution index\",\n",
    "    policy_index: \"Policy index\",\n",
    "    institution_index: \"Institutional index\"\n",
    "}\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 3. RUN MODELS\n",
    "# -------------------------------------------------\n",
    "table1_models = {}\n",
    "table1_summary_rows = []\n",
    "\n",
    "for dv in table1_dvs:\n",
    "    # Model 1: Other only\n",
    "    m1 = run_model(dv, include_emotion=False)\n",
    "    model_name_1 = f\"{dv}_Model1\"\n",
    "    table1_models[model_name_1] = m1\n",
    "\n",
    "    # Model 2: Other + Emotion\n",
    "    m2 = run_model(dv, include_emotion=True)\n",
    "    model_name_2 = f\"{dv}_Model2\"\n",
    "    table1_models[model_name_2] = m2\n",
    "\n",
    "    # Extract summary rows\n",
    "    for model_label, model in [(\"Model 1\", m1), (\"Model 2\", m2)]:\n",
    "        conf = model.conf_int()\n",
    "\n",
    "        row = {\n",
    "            \"DV\": dv_pretty[dv],\n",
    "            \"Model\": model_label,\n",
    "            \"N\": int(model.nobs),\n",
    "            \"R2\": model.rsquared,\n",
    "            \"Adj_R2\": model.rsquared_adj,\n",
    "        }\n",
    "\n",
    "        # other\n",
    "        if \"other\" in model.params.index:\n",
    "            row[\"Other_coef\"] = model.params[\"other\"]\n",
    "            row[\"Other_se\"] = model.bse[\"other\"]\n",
    "            row[\"Other_ci_low\"] = conf.loc[\"other\", 0]\n",
    "            row[\"Other_ci_high\"] = conf.loc[\"other\", 1]\n",
    "            row[\"Other_p\"] = model.pvalues[\"other\"]\n",
    "        else:\n",
    "            row[\"Other_coef\"] = None\n",
    "            row[\"Other_se\"] = None\n",
    "            row[\"Other_ci_low\"] = None\n",
    "            row[\"Other_ci_high\"] = None\n",
    "            row[\"Other_p\"] = None\n",
    "\n",
    "        # emotional burden\n",
    "        if emotion_index in model.params.index:\n",
    "            row[\"Emotion_coef\"] = model.params[emotion_index]\n",
    "            row[\"Emotion_se\"] = model.bse[emotion_index]\n",
    "            row[\"Emotion_ci_low\"] = conf.loc[emotion_index, 0]\n",
    "            row[\"Emotion_ci_high\"] = conf.loc[emotion_index, 1]\n",
    "            row[\"Emotion_p\"] = model.pvalues[emotion_index]\n",
    "        else:\n",
    "            row[\"Emotion_coef\"] = None\n",
    "            row[\"Emotion_se\"] = None\n",
    "            row[\"Emotion_ci_low\"] = None\n",
    "            row[\"Emotion_ci_high\"] = None\n",
    "            row[\"Emotion_p\"] = None\n",
    "\n",
    "        table1_summary_rows.append(row)\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 4. SUMMARY TABLE\n",
    "# -------------------------------------------------\n",
    "table1_df = pd.DataFrame(table1_summary_rows)\n",
    "\n",
    "# nicer rounding for display/save\n",
    "table1_display = table1_df.copy()\n",
    "num_cols = [\n",
    "    \"R2\", \"Adj_R2\",\n",
    "    \"Other_coef\", \"Other_se\", \"Other_ci_low\", \"Other_ci_high\", \"Other_p\",\n",
    "    \"Emotion_coef\", \"Emotion_se\", \"Emotion_ci_low\", \"Emotion_ci_high\", \"Emotion_p\"\n",
    "]\n",
    "for c in num_cols:\n",
    "    if c in table1_display.columns:\n",
    "        table1_display[c] = table1_display[c].round(3)\n",
    "\n",
    "# Save main table\n",
    "table1_display.to_excel(os.path.join(TAB_DIR, \"Table1_6models.xlsx\"), index=False)\n",
    "table1_display.to_csv(os.path.join(TAB_DIR, \"Table1_6models.csv\"), index=False, encoding=\"utf-8-sig\")\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 5. FULL RESULTS\n",
    "# -------------------------------------------------\n",
    "table1_full = save_full_results(\n",
    "    table1_models,\n",
    "    outfile_excel=os.path.join(RES_DIR, \"table1_6models_full_results.xlsx\"),\n",
    "    outfile_csv_prefix=os.path.join(RES_DIR, \"table1_6models_full_results\")\n",
    ")\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 6. SAVE TEXT SUMMARIES\n",
    "# -------------------------------------------------\n",
    "with open(os.path.join(TAB_DIR, \"Table1_6models_summaries.txt\"), \"w\", encoding=\"utf-8\") as f:\n",
    "    for name, model in table1_models.items():\n",
    "        f.write(f\"{'='*110}\\n\")\n",
    "        f.write(f\"{name}\\n\")\n",
    "        f.write(f\"{'='*110}\\n\\n\")\n",
    "        f.write(model.summary().as_text())\n",
    "        f.write(\"\\n\\n\\n\")\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 7. DISPLAY IN JUPYTER\n",
    "# -------------------------------------------------\n",
    "display(Markdown(\"## Table 1\"))\n",
    "display(Markdown(\"### Main results: 6 models\"))\n",
    "display(table1_display)\n",
    "\n",
    "display(Markdown(\"### Full regression results\"))\n",
    "display(table1_full)\n",
    "\n",
    "display(Markdown(\"### Model summaries\"))\n",
    "\n",
    "for name, model in table1_models.items():\n",
    "    display(Markdown(f\"#### {name}\"))\n",
    "    print(model.summary())\n",
    "    print(\"\\n\" + \"=\"*120 + \"\\n\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4938d589-4ecd-4834-bd32-d173c328824a",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "3409392c-cc60-49b1-9da6-5b8476a6fe42",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "id": "3bdf5378-ea40-4835-8b72-32ddbec1b196",
   "metadata": {},
   "source": [
    "## Figure 3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "cbf11161-3ab2-4845-a237-1d2cb4097135",
   "metadata": {},
   "outputs": [],
   "source": [
    "# -------------------------------------------------\n",
    "# 1. OUTCOME SETUP\n",
    "# -------------------------------------------------\n",
    "appendix_outcomes = [\n",
    "    attribution_index,\n",
    "    policy_index,\n",
    "    institution_index\n",
    "]\n",
    "\n",
    "appendix_titles = {\n",
    "    attribution_index: \"(a) Predicted attribution\",\n",
    "    policy_index: \"(b) Predicted policy support\",\n",
    "    institution_index: \"(c) Predicted institutional evaluation\"\n",
    "}\n",
    "\n",
    "appendix_ylabels = {\n",
    "    attribution_index: \"Predicted attribution index\",\n",
    "    policy_index: \"Predicted policy index\",\n",
    "    institution_index: \"Predicted institutional index\"\n",
    "}\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 2. RUN INTERACTION MODELS\n",
    "# -------------------------------------------------\n",
    "appendix_models = {}\n",
    "appendix_predicted_all = []\n",
    "appendix_key_rows = []\n",
    "\n",
    "for dv in appendix_outcomes:\n",
    "    m = run_model(dv, include_asset_interaction=True)\n",
    "    model_name = f\"Appendix_{dv}_interaction\"\n",
    "    appendix_models[model_name] = m\n",
    "\n",
    "    # -----------------------------\n",
    "    # Save key terms for compact table\n",
    "    # -----------------------------\n",
    "    conf = m.conf_int()\n",
    "    for term in [\"other\", asset_index, f\"other:{asset_index}\"]:\n",
    "        if term in m.params.index:\n",
    "            appendix_key_rows.append({\n",
    "                \"model\": model_name,\n",
    "                \"dv\": dv,\n",
    "                \"term\": term,\n",
    "                \"coef\": m.params[term],\n",
    "                \"se\": m.bse[term],\n",
    "                \"ci_low\": conf.loc[term, 0],\n",
    "                \"ci_high\": conf.loc[term, 1],\n",
    "                \"p_value\": m.pvalues[term],\n",
    "                \"n\": int(m.nobs),\n",
    "                \"r2\": m.rsquared,\n",
    "                \"adj_r2\": m.rsquared_adj\n",
    "            })\n",
    "\n",
    "    # -----------------------------\n",
    "    # prediction grid\n",
    "    # -----------------------------\n",
    "    asset_grid = np.linspace(df[asset_index].quantile(0.02), df[asset_index].quantile(0.98), 120)\n",
    "\n",
    "    pred_base = {}\n",
    "    for var in continuous_controls:\n",
    "        pred_base[var] = df[var].mean()\n",
    "\n",
    "    for var in categorical_controls:\n",
    "        pred_base[var] = df[var].mode(dropna=True)[0]\n",
    "\n",
    "    pred_rows = []\n",
    "    for other_val in [0, 1]:\n",
    "        for a in asset_grid:\n",
    "            row = pred_base.copy()\n",
    "            row[\"other\"] = other_val\n",
    "            row[asset_index] = a\n",
    "            pred_rows.append(row)\n",
    "\n",
    "    pred_df_tmp = pd.DataFrame(pred_rows)\n",
    "    pred_tmp = m.get_prediction(pred_df_tmp).summary_frame(alpha=0.05)\n",
    "\n",
    "    pred_df_tmp[\"pred\"] = pred_tmp[\"mean\"].values\n",
    "    pred_df_tmp[\"ci_low\"] = pred_tmp[\"mean_ci_lower\"].values\n",
    "    pred_df_tmp[\"ci_high\"] = pred_tmp[\"mean_ci_upper\"].values\n",
    "    pred_df_tmp[\"condition\"] = np.where(pred_df_tmp[\"other\"] == 0, \"Self\", \"Other\")\n",
    "    pred_df_tmp[\"dv\"] = dv\n",
    "    pred_df_tmp[\"dv_label\"] = item_labels[dv]\n",
    "    pred_df_tmp[\"model\"] = model_name\n",
    "\n",
    "    appendix_predicted_all.append(pred_df_tmp)\n",
    "\n",
    "appendix_predicted_df = pd.concat(appendix_predicted_all, ignore_index=True)\n",
    "appendix_key_df = pd.DataFrame(appendix_key_rows)\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 3. SAVE FULL RESULTS\n",
    "# -------------------------------------------------\n",
    "appendix_reg_full = save_full_results(\n",
    "    appendix_models,\n",
    "    outfile_excel=os.path.join(RES_DIR, \"appendix_figureA1_interaction_full_results.xlsx\"),\n",
    "    outfile_csv_prefix=os.path.join(RES_DIR, \"appendix_figureA1_interaction_full_results\")\n",
    ")\n",
    "\n",
    "appendix_predicted_df.to_excel(\n",
    "    os.path.join(RES_DIR, \"appendix_figureA1_predicted_values.xlsx\"),\n",
    "    index=False\n",
    ")\n",
    "appendix_predicted_df.to_csv(\n",
    "    os.path.join(RES_DIR, \"appendix_figureA1_predicted_values.csv\"),\n",
    "    index=False,\n",
    "    encoding=\"utf-8-sig\"\n",
    ")\n",
    "\n",
    "# compact key-terms results\n",
    "appendix_key_df.to_excel(\n",
    "    os.path.join(RES_DIR, \"appendix_figureA1_key_terms.xlsx\"),\n",
    "    index=False\n",
    ")\n",
    "appendix_key_df.to_csv(\n",
    "    os.path.join(RES_DIR, \"appendix_figureA1_key_terms.csv\"),\n",
    "    index=False,\n",
    "    encoding=\"utf-8-sig\"\n",
    ")\n",
    "\n",
    "# save text summaries\n",
    "with open(os.path.join(RES_DIR, \"appendix_figureA1_model_summaries.txt\"), \"w\", encoding=\"utf-8\") as f:\n",
    "    for name, model in appendix_models.items():\n",
    "        f.write(\"=\" * 120 + \"\\n\")\n",
    "        f.write(f\"{name}\\n\")\n",
    "        f.write(\"=\" * 120 + \"\\n\\n\")\n",
    "        f.write(model.summary().as_text())\n",
    "        f.write(\"\\n\\n\\n\")\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 4. DRAW FIGURE\n",
    "# -------------------------------------------------\n",
    "fig, axes = plt.subplots(\n",
    "    nrows=1, ncols=3,\n",
    "    figsize=(16.8, 6.5),\n",
    "    sharex=True,\n",
    "    sharey=False\n",
    ")\n",
    "\n",
    "for ax, dv in zip(axes, appendix_outcomes):\n",
    "    temp_dv = appendix_predicted_df[appendix_predicted_df[\"dv\"] == dv].copy()\n",
    "\n",
    "    for cond, ls, marker in [(\"Self\", \"-\", \"o\"), (\"Other\", \"--\", \"s\")]:\n",
    "        temp = temp_dv[temp_dv[\"condition\"] == cond].copy()\n",
    "\n",
    "        ax.plot(\n",
    "            temp[asset_index],\n",
    "            temp[\"pred\"],\n",
    "            linestyle=ls,\n",
    "            linewidth=2.0,\n",
    "            marker=None,\n",
    "            color=\"black\",\n",
    "            label=cond,\n",
    "            zorder=3\n",
    "        )\n",
    "\n",
    "        ax.fill_between(\n",
    "            temp[asset_index],\n",
    "            temp[\"ci_low\"],\n",
    "            temp[\"ci_high\"],\n",
    "            alpha=0.12,\n",
    "            color=\"black\",\n",
    "            zorder=2\n",
    "        )\n",
    "\n",
    "        # anchor markers\n",
    "        anchor_idx = np.linspace(0, len(temp) - 1, 6, dtype=int)\n",
    "        temp_anchor = temp.iloc[anchor_idx]\n",
    "        ax.plot(\n",
    "            temp_anchor[asset_index],\n",
    "            temp_anchor[\"pred\"],\n",
    "            linestyle=\"None\",\n",
    "            marker=marker,\n",
    "            markersize=6.3,\n",
    "            color=\"black\",\n",
    "            zorder=4\n",
    "        )\n",
    "\n",
    "        # annotate first/last points\n",
    "        annotate_endpoints(\n",
    "            ax,\n",
    "            temp[asset_index].values,\n",
    "            temp[\"pred\"].values,\n",
    "            condition=cond,\n",
    "            label_prefix=\"\",\n",
    "            fontsize=14.0,\n",
    "            y_offset=0.05\n",
    "        )\n",
    "\n",
    "    ax.set_title(appendix_titles[dv], loc=\"center\", pad=TITLE_PAD, fontsize=18)\n",
    "    ax.set_xlabel(\"Asset index\")\n",
    "    ax.set_ylabel(appendix_ylabels[dv])\n",
    "\n",
    "    ax.tick_params(axis='both', which='major', length=5, width=1.0, direction='out')\n",
    "    ax.tick_params(axis='both', which='minor', length=3, width=0.8, direction='out')\n",
    "    ax.spines[\"left\"].set_linewidth(1.0)\n",
    "    ax.spines[\"bottom\"].set_linewidth(1.0)\n",
    "\n",
    "    ax.legend(\n",
    "        loc=\"best\",\n",
    "        frameon=False,\n",
    "        handlelength=1.9,\n",
    "        borderpad=0.15,\n",
    "        labelspacing=0.28,\n",
    "        handletextpad=0.4,\n",
    "        fontsize=14\n",
    "    )\n",
    "\n",
    "plt.tight_layout()\n",
    "fig.savefig(os.path.join(FIG_DIR, \"Appendix_FigureA1_three_index_predicted.png\"), bbox_inches=\"tight\")\n",
    "fig.savefig(os.path.join(FIG_DIR, \"Appendix_FigureA1_three_index_predicted.pdf\"), bbox_inches=\"tight\")\n",
    "plt.show()\n",
    "\n",
    "# -------------------------------------------------\n",
    "# 5. DISPLAY FULL RESULTS IN JUPYTER\n",
    "# -------------------------------------------------\n",
    "display(Markdown(\"## Appendix Figure A1\"))\n",
    "\n",
    "display(Markdown(\"### Predicted values for attribution, policy, and institutional indices\"))\n",
    "display(appendix_predicted_df)\n",
    "\n",
    "display(Markdown(\"### Key interaction terms\"))\n",
    "display(appendix_key_df)\n",
    "\n",
    "display(Markdown(\"### Full regression results for Appendix Figure A1 interaction models\"))\n",
    "display(appendix_reg_full)\n",
    "\n",
    "display(Markdown(\"### Regression summaries\"))\n",
    "for name, model in appendix_models.items():\n",
    "    display(Markdown(f\"#### {name}\"))\n",
    "    print(model.summary())\n",
    "    print(\"\\n\" + \"=\" * 120 + \"\\n\")\n",
    "\n",
    "print(\"Saved files:\")\n",
    "print(os.path.join(FIG_DIR, \"Appendix_FigureA1_three_index_predicted.png\"))\n",
    "print(os.path.join(FIG_DIR, \"Appendix_FigureA1_three_index_predicted.pdf\"))\n",
    "print(os.path.join(RES_DIR, \"appendix_figureA1_predicted_values.xlsx\"))\n",
    "print(os.path.join(RES_DIR, \"appendix_figureA1_predicted_values.csv\"))\n",
    "print(os.path.join(RES_DIR, \"appendix_figureA1_key_terms.xlsx\"))\n",
    "print(os.path.join(RES_DIR, \"appendix_figureA1_key_terms.csv\"))\n",
    "print(os.path.join(RES_DIR, \"appendix_figureA1_interaction_full_results.xlsx\"))\n",
    "print(os.path.join(RES_DIR, \"appendix_figureA1_interaction_full_results.csv\"))\n",
    "print(os.path.join(RES_DIR, \"appendix_figureA1_model_summaries.txt\"))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "aef09152-6a05-4054-97b5-427d744ff37c",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1478ae3c-2207-4576-97a9-694dc6a1797f",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "78d9f285-8f50-4cff-beaf-1f4c3e21851d",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a67a49af-8d10-48a4-869f-a7898b7a7e02",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "327c7f17-6eb4-44db-8d20-5b08d2682528",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "9f65b8a0-a593-4150-9d40-906b6d80f4e1",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5df556f3-ac28-490e-8fba-1e26a8193207",
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f47508b7-b15f-4d34-8220-690a6d11d870",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.12.8"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
